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To constrain theories for the dynamical evolution of global ice mass through the 
late Neogene, it is important to determine whether major changes in the record were 
gradual or rapid. Of particular interest is the evolution of the near 100-ky ice age 
cycle in the middle Pleistocene. We have applied a new technique based on multiple 
taper spectrum analysis which allows us to model the time evolution of quasi- periodic 
signals [Park and Maasch, 1992]. This technique uses both phase and amplitude 
information, and enables us to address the question of abrupt versus gradual onset 
of the 100-ky periodicity in the middle Pleistocene. 

The variation of climate proxy variables at a given frequency f 0 (with associated 
period T 0 = 1 //„) can be parameterized by 3£{A(/.)e -2ir,/ot }, a sinusoid with a slowly- 
varying amplitude A (t). The function T(t) is complex-valued, allowing slow variations 
in phase as well as magnitude. We define A(t) as the ‘envelope 1 of the signal at ‘carrier 
frequency 1 f Q . If A(f) = A 0 is a constant, the signal is termed ‘periodic 1 or ‘a phase- 
coherent sinusoid. 1 If A(t) varies, the signal is termed ‘quasi-pcriodic. 1 All methods 
of determining the envelope of a quasi-periodic signal have shortcomings. Bandpass 
filters, such as those used by Ruddiman et al [1989], cannot model discontinuities in 
the envelope of a quasi-periodic signal, and thus, are of little help in discriminating 
between a sudden and a gradual onset lor the 100-ky ice age cycle. Moreover, a 
narrow bandpass in the frequency domain requires a long time-domain filter, so that 
a significant number of data points at the ends of the series must be discarded. 
Similarly, the choice of length for the overlapping data segments in the Ruddiman et 
al [1986b] analysis involves a tradeoff between frequency and time resolution. 

We apply instead a more flexible approach based on multiple-taper spectrum 
analysis [Thomson 1982; 1990; Park et al. 1987, Berger et al 1991]. The multitaper 
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approach allows the analyst to model the envelope function A(t) using the tools of lin- 
ear inverse theory, solving for the envelope function that fits the time series data while 
optimizing some property of the envelope. Inversion algorithms can be derived that 
find the ‘smallest’ and ‘smoothest’ envelopes that fit spectral ‘data’ from a given time 
series [Park 1990; Park 1992], Special cases of these the algorithms can be derived 
to model sudden changes in the envelope. This allows us to examine the abruptness 
of the onset of the 100-ky periodicity as well as the evolution of the obliquity and 
precession signals. The estimation procedure relates spectral estimates at / 0 , using a 
set of orthogonal Slepian data tapers, as linear functionals of an unknown envelope 
/1(f), so that the envelope is constructed from a linear combination of the Slepian 
tapers. Since the Slepian tapers are optimally bandlimited, multitaper envelope es- 
timation can be thought oi as an extremely sharp narrow-band filter valid for the 
entire duration of the series. 

The shortcomings of the algorithm we use are similar to those of other linear 
inverse problems. For instance, the multitaper envelope estimation procedure can 
model envelope discontinuities, but the technique cannot by itself discriminate be- 
tween continuous and discontinuous models for the 100-ky cycle, since envelopes of 
both types can be constructed that fit the spectral ‘data’ exactly. The analyst must 
use other, perhaps subjective, criteria for choosing among models that fit the data. 
I his nonuniqueness is a common problem in linear inverse theory; in principle an 
infinite number of envelopes can fit a finite number of spectral data exactly. In 
most cases, the a priori constraint that the signal at a carrier frequency f Q have a 
‘slowly- varying’ envelope, except at one or more specified time points, reduces greatly 
the numbei of acceptable models for the envelope, and motivates the solution of a 
‘smoothest’ envelope that fits the spectral data. 

We investigated the time-evolution of the 100-ky cycle in 6 m O data, thought to 
reflect global ice-volume variations. We analyzed 6 18 0 data from DSDP Site 607 
and ODP Site 677, from which three long (> 2.6 My) time series have been published 
[Shackleton and Hall, 1989; Ruddiman et al 1986ab; 1989; Raymoet al 1989]. We find 
evidence for a coherent <5 18 0 signal for both cores in the eccentricity and obliquity 
frequency bands, consistent with variations in global ice volume as the causative 
factor. However, the nature of the earth system response to orbital insolation cycles 
depends on the time scale adopted in the spectral analysis. If the Ruddiman/Raymo 
time scale for DSDP Site 607 is accepted, the <*> 18 0 obliquity cycle has enhanced 
amplitude between 1.0 and 1.5 Ma, relative to the late Pleistocene (t < 1.0 Ma), and 
a nonlinear response of the earth system to orbital obliquity is inferred (Figure 1). 
If the Shackleton et al [1990] time scale for ODP Site 677 is accepted, the amplitude 
match between the S 18 0 obliquity cycle and the 65°N insolation derived by Berger 
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and Loutre [1991] from the recent astronomical solution of Laskar [1988; 1990] is very 
good for times t< 2.3 Ma (Figure 2). (The phase match for obliquity is enhanced by 
the fact that the data series were tuned to the astronomical time series, so that the 
observed phase agreement is not surprising.) The veracity of the ODP 677 time scale 
therefore correlates with a linear earth- system response to orbital obliquity. 

Based on our analysis of data from these two deep-sea cores, we do not find 
compelling evidence for an abrupt change in the 100-ky <5 18 0 signal. Rather, envelope 
inversions in the eccentricity band suggest that the 100-ky 6 ls O cycle is phase-locked 
with the 124-ky eccentricity cycle some 300-400 ky prior to its late-Pleistocene growth 
in amplitude and phase-lock with the 95-ky eccentricity cycle (Figure 3). An abrupt 
change in the 95-ky envelope near 0.85 Ma is consistent with DSDP 607 data on the 
Ruddiman/Raymo time scale, but such a transition would occur while leaving the 
124-ky envelope largely unchanged. If the Shackleton et al [1990] time scale for ODP 
677 is accepted, our three 8 18 0 records are consistent with a low-amplitude 100-ky 
cycle between 1.2 and 2.6 Ma, whose local period of oscillation alternates between the 
two major eccentricity periods at 95 and 124 ky. The phase of these 100-ky oscillations 
prior to 1.2 Ma is related to the phase of the astronomical eccentricity cycles where 
the cycles have significant amplitude. The match between the precessional envelopes 
of Sites 607 and 677 is poor, when both are expressed in terms of the Shackleton et al 
[1990] time scale. Climate simulation studies suggest that cyclic salinity changes in 
equatorial surface waters are a plausible contributor to the ODP 677 <5 18 0 data in the 
precessional band, and could explain this discrepancy. However, our time- rescaling of 
data from DSDP 607, using visual isotope-stage matching, may not be precise enough 
to transform the short- period precession cycles with sufficient accuracy. 

Further study of more paleoclimate records will be necessary in order to address 
more fully the time-evolution of the 100 kyr cycle. For instance, the correlation of 
£ 18 0 signals from Sites 677 and 607 is intriguing, but comprehensive tests for the 
global coherence of <5 18 0 signals should use data from more than two sites. The 
phase information contained in the paleoclimate records may reaffirm the notion that 
external earth-orbital forcing could be the pacemaker of the ice age cycles. How- 
ever, the mid-Pleistocene amplitude increase of the quasi-periodic 100-ky <5 18 0 signal 
awaits a complete explanation, which appears to depend on factors other than orbital 
insolation changes. 

Future Tasks: 

1) Collection of climate proxy data series, with good time control, from both 
deep-sea, lake-bed, and land-based sedimentary sequences. Much data probably lies 
dormant in the older cores stored by the Ocean Drilling Program, but only the older 
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cores often suffer from coring gaps that impede detailed time series analysis. Newer 
core data are spliced from parallel-paired drill cores at a single location. 

2) careful cross-spectral analysis of different data series to investigate how the 
earth climate system, as a whole, responded to orbital insolation cycles over the past 
2-3 My, the period of Earth history most relevant to current global change problems. 
The techniques developed for the above project can be extended to cross-series studies. 
The buzzword for (his kind of study is ’global teleconnections, ’ the manner in which 
climate at different parts of the globe interacts. The response of climate to modest 
changes in its boundary conditions can (in principle) be tracked by its response over 
individual Milankovitch cycles. 

3) similar analysis of data series from earlier periods of Earth history e.g. the 
Cretaceous and the Eocene. Data is accumulating for a significant response to orbital 
cycles, but a global synthesis has only been attempted so far with numerical climate 
models. The Earth’s climate was warmer in these two periods, but our understanding 
of how the Earth sustained such temperature is, at best, incomplete. Many studies 
point towards higher carbon dioxide levels as a causative factor, and the ubiquitous 
Milankovich-driven limestone/black-shale sequences suggest big changes in ocean cir- 
culation. Are the Cretaceous and Eocene global climates a picture of what awaits us 
in a greenhouse future? 

4) Despite a decade of global-wanning predictions based on numerical climate 
models, it must be admitted that these models (global circulation, energy balance, 
etc) can represent the earth’s climate dynamics only in a crude manner. The features 
of the climate system that are critical to its longer-term variation (decades and cen- 
turies) may not be apparent with the current generation of numerical climate models. 
Improving the global circulation models is a high priority. This includes atmospheric, 
ocean and coupled atmosphere-ocean models. 

5) Collection of a long-term global climate database to ’validate’ the output of 
global circulation models. Such data are essential to check if the numerical climate 
models are working. Comparisons of atmospheric GCM results with satellite data 
have begun, but I am not aware of a global data- validation of oceanic GCMs. 
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Figure Captions 


Figure 1. Envelope inversions for the time-evolution of the 41-kv signal in the 
benthic 6 l8 0 series from DSDP Site 607 and ODP Site 677, plotted against similar 
analyses lor the astronomical insolation series derived from Berger f 19781 (‘65° old’) 
and Berger and Loutre [1991] (‘65°N new’). Seven 47r-prolate Slepian eigentapers 
aie used to constrain the estimates, using ’smoothness’ penalty function (10). The 
amplitude units in this and succeeding plots are parts-per-thousand d) 18 0 Note the 
phase * correlation (with a constant shift) between the older astronomical series and 
data from DSDP Site 607, and the correlation between the newer astronomical series 
and data from ODP Site 677. This arises from the orbital tuning performed on the 
data series. Orbital tuning does not enforce correlations between envelope amplitudes 
as is evident from the upper panel. ’ 


Figure 2. Similar to Figure 1, but with the envelope of the Berger [1978] se- 
ries omitted and benthic <5 18 0 data from DSDP Site 607 expressed in terms of the 
S hacklet on et al [1990] time scale. Seven 47r-prolate Slepian eigentapers are used to 
constrain the estimates, using a ’smoothness’ penalty function. Note the phase corre- 
lation (with a constant shift) between the data series. Note the improved correlations 
between envelope amplitudes in the upper panel. 


Figure 3. 1 est for the abrupt onset of the 100-ky signal. Double-line envelope 
inversion for the benthic/ 18 0 series from DSDP Site 607 on its published time scale 
at the two major eccentricity quasi-periods of 95 and 124 ky. Seven 47r-prolate Slepian 
eigentapers are used to constrain the envelopes, which are constrained to be smooth 
everywhere except at t 0 =0.85 Ma. 
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Envelopes for Obliquity Cycle 
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